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The calculation of a highly accurate cosmological recombination history has been the object of 
particular attention recently, as it constitutes the major theoretical uncertainty when predicting 
the angular power spectrum of Cosmic Microwave Background anisotropies. Lyman transitions, in 
particular the Lyman-a line, have long been recognized as one of the bottlenecks of recombination, 
due to their very low escape probabilities. The Sobolev approximation does not describe radiative 
transfer in the vicinity of Lyman lines to a sufficient degree of accuracy, and several corrections 
have already been computed in other works. In this paper, the impact of some previously ignored 
radiative transfer effects is calculated. First, the effect of Thomson scattering in the vicinity of the 
Lyman-a line is evaluated, using a full redistribution kernel incorporated into a radiative transfer 
code. The effect of feedback of distortions generated by the optically thick deuterium Lyman-a 
line blueward of the hydrogen line is investigated with an analytic approximation. It is shown that 
both effects are negligible during cosmological hydrogen recombination. Secondly, the importance of 
high-lying, non overlapping Lyman transitions is assessed. It is shown that escape from lines above 
Ly7 and frequency diffusion in Ly/3 and higher lines can be neglected without loss of accuracy. 
Thirdly, a formalism generalizing the Sobolev approximation is developed to account for the overlap 
of the high-lying Lyman lines, which is shown to lead to negligible changes to the recombination 
history. Finally, the possibility of a cosmological hydrogen recombination maser is investigated. It 
is shown that there is no such maser in the purely radiative treatment presented here. 



I. INTRODUCTION 

The first measurements of the Cosmic Microwave 
Background (CMB) spectrum pQ and temperature 
anisotropies [2] changed cosmology from a qualitative to 
a robust and predictive science. Since then our picture of 
the Universe has become more and more accurate. Ob- 
servations of high-redshift type la supernovae [21 H] have 
made it clear that nearly three fourths of the energy bud- 
get of our Universe is a non-clustering "dark energy" fluid 
with a negative pressure. In the last decade, the measure- 
ments of the temperature and polarization anisotropies in 
the CMB by the Wilkinson Microwave Anisotropy Probe 
(WMAP) [5] have confirmed this picture and propelled 
cosmology into the era of high precision. Combined with 
other CMB measurements (e.g. BOOMERANG [5J, CBI 
[7j, ACBAR [H], QUaD [§]) and large-scale structure sur- 
veys (2dF [10!, SDSS 03]), WMAP results have firmly 
established the ACDM model as the standard picture of 
our Universe. 

What has also emerged from this high-precision data is 
our ignorance of the large majority of the constituents of 
the Universe. Only ~ 4% of our Universe is in the form of 
known matter (most of which is not luminous), the rest 
is in the form of an unknown clustering "dark matter" 
(~ 23%) or the even more disconcerting "dark energy" 
(~ 73%). In addition, it is now widely believed that 
the Universe underwent an inflationary phase early on 
that sourced the nearly scale-invariant primordial density 
perturbations, which led to the large-scale structure we 
observe today. Inflation requires non-standard physics, 
and at present there is no consensus on the mechanism 
that made the Universe inflate, and only few constraints 



on the numerous inflationary models are available from 
observations. 

The Planck satellite, launched in May 2009, will mea- 
sure the power spectrum of temperature anisotropies in 
the CMB, Cj T , with a sub-percent accuracy, up to mul- 
tipole moments t ~ 2500 |12j . It will also measure 
the power spectrum of E-mode polarization anisotropies 
up to £ ~ 1500. With this unprecedented ultra-high- 
precision data, cosmologists will be in a position to in- 
fer cosmological parameters accurate to the sub-percent 
level. The high resolution of Planck observations will 
provide a lever arm to precisely measure the spectral in- 
dex of scalar density perturbations n§ and their running 
as, therefore usefully constraining models of inflation. 
The polarization data will help break degeneracies of cos- 
mological parameters with the optical depth to the sur- 
face of last scattering r, giving us a better handle on 
the epoch of reionization. This wealth of upcoming high- 
precision data from Planck, as well as that from ongo- 
ing experiments (ACT [T3], SPT [13] ) or possible future 
space-based polarization missions (CMBPol [H]), can be 
fully exploited only if our theoretical predictions of CMB 
anisotropies are at least as accurate as the data. 

The physics of CMB anisotropy generation is now well 
understood, and public Boltzmann codes are available 
(CMBFast [16], CAMB [T7], CMBEasy [H]), which 
evolve the linear equations of matter and radiation per- 
turbations and output highly accurate CMB temperature 
and polarization angular power spectra, for a given ion- 
ization history 19]. The dominant source of systematic 
uncertainty in the predicted C^s is the recombination his- 
tory [20]. Not only the peak and width of the visibility 
function are important, but the precise shape of its tails 
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is also critical at the sub-percent level of accuracy, in 
particular for the Silk damping tail |21j of the anisotropy 
power spectrum. This has motivated Seager et al. [2"2"ll2"3"] 
to revise the seminal work of Peebles [21] and Zeldovich 
et al. 25 and extend their three-level atom model to 
a multi-level atom (MLA) calculation. Their recombina- 
tion code, RecFast, is accurate to the percent level, and 
is a part of the Boltzmann codes routinely used for cur- 
rent day CMB data analysis. While sufficiently accurate 
for WMAP data, RecFast does not satisfy the level of 
accuracy required by Planck [25] [27] . 

In the last few years, significant work has been de- 
voted to further understanding the rich physics of cos- 
mological recombination. Helium recombination is now 
understood to a sufficient level of accuracy for Planck 
[25H55] . Cosmological hydrogen recombination demands 
a much higher level of accuracy than helium [36] , and its 
precise description is still the subject of ongoing efforts. 
Previous work on hydrogen recombination can be cast 
into two categories. 

On the one hand, accurate recombination histories 
need to account for as large a number as possible of ex- 
cited states of hydrogen. This is particularly important 
at late times, z < 800 — 900, when the free electron abun- 
dance becomes very low and the slow recombinations to 
the excited states become the "bottleneck" of the recom- 
bination process. In these conditions, it is important to 
precisely account for all possible recombination pathways 
by including a large number of excited states in MLA cal- 
culations. Since the recombination rates strongly depend 
on the angular momentum quantum number I, an accu- 
rate code must resolve the angular momentum substates 
[57t [38] and lift the statistical equilibrium assumption 
previously made. The standard MLA approach requires 
solving for the population of all the excited states ac- 
counted for, which is computationally expensive and has 
limited recent high-rt computations 39, 40J to only a few 
points in parameter space. Recently, a new effective MLA 
(EMLA) method has been introduced [UJ, which makes 
it possible to account for a very large number of excited 
states, while preserving the computational efficiency of a 
simple few-level atom. The EMLA approach consists in 
factoring the effect of the "interior" excited states (states 
which are not connected to the ground state) into effec- 
tive recombination and photoionization coefficients and 
bound bound transition rates for the small number of "in- 
terface" states radiatively connected to the ground state, 
i.e. 2s, 2p and the low-lying p states. 

The second category of studies has concentrated on 
the transitions from the "interface" states to the ground 
state, of particular importance at early times when the 
overall recombination rate is controlled by the so called 
"n=2 bottleneck" . Hydrogen atoms in the n~2 shell can 
reach the ground state either through emissions of sin- 
gle photons from the 2p state, which slowly escape the 
optically thick Lyman-a transition through cosmological 
redshifting, or from the 2s state, through forbidden two- 
photon decays. Previous studies have examined a series 



of effects that may affect these transition rates, including 
but not limited to stimulated 2s — > Is two-photon decays 
and non-thermal Is — > 2s two-photon absorptions |42L 
HI], feedback between neighboring lines of the optically 
thick Lyman series [HI US] , time-dependent effects in the 
Lyman-a line |44[ I47j , two-photon decays from higher ex- 
cited states [1H H5H50] and diffusion in the Lyman-a line 

IEH53]- 

The purpose of this work is mainly to assess the im- 
portance of a few radiative transfer effects in the Lyman 
lines, that have not been investigated yet, or not in suffi- 
cient detail. In Section [H] we review the theory of radia- 
tive transfer in the vicinity of a resonant line in an ex- 
panding Universe. We then turn to the Lyman-a line, for 
which we study the effect of Thomson scattering (Section 



III A ) and the interaction with the Lya line of deuterium 
(Section III B ) . In Section IV we quantify the impor- 



tance of the high-lying, non-overlapping Lyman transi- 
tions. Section [V] is devoted to the overlapping high- lying 
Lyman lines. We also explore the possibility of a cosmo- 
logical recombination maser in Sec tion|VI| We summarize 
and discuss our results in Section |VII[ In this paper we 
will not repeat the standard MLA nor the EMLA formal- 
ism, and refer the reader to Refs. [22[ HI] and [IT] , 
respectively, for a detailed description. 



II. RADIATIVE TRANSFER IN THE LYMAN 
LINES 

A. Basic notation 

In this section, we present the basic quantities and 
notation used throughout this paper. We summarize our 
notation in Table [TJ 

The photon occupation number at frequency v is de- 
noted /„. In the case of a black body spectrum with 
temperature T, /„ = (e hv / kT — i) . We will also make 
use of the number of photons per hydrogen atom per unit 
frequency, 



^r 8 ™ 2 , 
c 6 Nh 



(1) 



where Nr is the number density of hydrogen atoms. The 
population of a species X relative to the total abundance 
of hydrogen is denoted ix = -Vx/^Vh- The fractional 
abundance of hydrogen in the state [n, I] is denoted x n i- 
For the low I states, we use the spectroscopic notation 
s,p, d, so the ground state is denoted Is and the [n, I = 
1] states are denoted np. 

This work will be concerned primarily with the np — > 
Is transitions, which will be referred to as the Lyman- 
n (or Ly-n) transitions. The Ly-2 transition therefore 
designates, in that convention, the Lyman-a (Lya) tran- 
sition. Wc denote the resonant Ly-n transition frequency 



1 - 



1 



"Lya = 1 - 



(2) 
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TABLE I. Notation used in this paper. Units of "1" means the quantity is dimensionless. 



Symbol Units 



Description 



Equation 
Eq. Q 

a n = r np /(47r^ ri A H ) 

/cq = X n v I (3 ^ls ) 

Eq. 1 15 I 



Anp, Is 
<X„ 

u 

rn 
.1 eq 

rrt 

/(em) 

H 

N u 
N x 
M v 

p: b 

■^n 1 l—>np 
Rn{v, V') 

R T {u v') 
11(5) 

s 

s 

T T 

w„ 
vv 
ir 

■''D 

X (: 
X 



S 

1 
1 
1 
1 

s" 1 
cm -3 
cm" 3 
Hz" 1 

1 

1 
s- 1 
Hz" 2 
Hz- 1 

1 

Hz 
Hz 
1 

K 
Hz 
Hz 

1 

1 

1 

1 

1 



X n 


1 


y 


1 


;/d 


1 


Trap 


s" 1 


A x 


1 




Hz 




1 


'/e 


Hz" 1 



Hyn (^H,^d) 


Hz 


Vc 


Hz 


Vn 


Hz 


0"t 


cm 2 


(TO 


cm 2 


T n 


1 




1 


T c 


1 




Hz" 1 


cj>v(x; a) 


1 


<t>n{x) 


1 


<t>o(x) 


1 


* (H) (y) 


1 


*(H+D) (y) 


1 




1 




1 



Ju 



C 3JV H 
Psc — ^np, Is f~^np 
Pab = l-?& 



Einstein A-coemcient for the spontaneous decay np — > Is 
dimensionless Voigt parameter for the Ly-n transition 
photon occupation number at frequency v 
equilibrium value of f v near Ly-n 
equilibrium value of /„ near Ly-n, with true absorption and emision only 
Hubble expansion rate 
total number density of hydrogen nuclei 
number density of species X 
number of photons per hydrogen atom per unit frequency M v 
fraction of photon absorptions in Ly-n resulting in a scattering 
fraction of true photon absorptions in Ly-n 
radiative transition rate n'l — > np per hydrogen atom in the [n', I] state 

resonant scattering redistribution kernel in Ly-n Eq. ( 28 1 

Thomson scattering kernel Eq. (51 1 

dimensionless Thomson scattering kernel Eq. (56 1 

width over which resonant scattering is effective in Ly-n Eq. ( 36 1 

width over which resonant scattering is effective in Lya S = S2 

dimensionless width over which resonant scattering is effective in Lya S = (hS/(kT T )) 3 

matter and radiation temperatures 
width over which the Ly-n wings are optically thick for true absorption 
width over which the Lya wings are optically thick for true absorption 
dimensionless width over which Lya is optically thick for true absorption 
abundance of deuterium relative to total hydrogen nuclei 
fraction of hydrogen in the state [n, I] 
abundance of free electrons relative to total hydrogen nuclei 
detuning from line center (when considering a single line), 
or from v c (when considering line overlap), in Doppler units 
detuning between Ly-n and Ly-continuum frequencies in Doppler widths 
same as x OR detuning from the Lya frequency in units of kT T /h 
detuning of the H and D Lya frequencies in units of kT T /h 
width (or inverse lifetime) of the np state 
dimensionless Doppler width for a scattering species of mass mx 
rms frequency shift during an electron scattering event near Lya 

relative amplitude of the distortion due to deuterium at D Lya Eq. (74 1 

differential optical depth for Thomson scattering near Lyman-a Eq. ( 50 1 

Lyman-a resonant frequency (in hydrogen or deuterium specifically) 

Lyman-limit frequency v c = ^VL.yc 

Lyman-n resonant frequency v n — (1 — n~' 

Thomson cross section 
photoionization cross section from the ground state, at threshold Eq. ( 99 1 

Sobolev optical depth in the Lyman-n resonance line Eq. (jfsj) 

Sobolev optical depth in the Lyman-a line r Lya = t~2 

Lyman continuum optical depth per Doppler width near threshold Eq. ( 111 I 

line profile for the Ly-n transition 
dimensionless Voigt profile, with Voigt parameter a 
dimensionless Doppler profile centered at the Ly-n frequency 
non-normalized dimensionless profile for continnum absorption near threshold Eq. ( 103 1 

dimensionless photon occupation number near Lya whith hydrogen only Eqs. (66 1, (' 

dimensionless photon occupation number near Lya whith deuterium Eqs. (66 1, (711 

dimensionless equilibrium photon occupation number near D Lya Eq^J 72 1 

rescaled difference in photon occupation numbers with or without deuterium Eqs. (751, (761 



Eq. |35l 

W = W 2 
W = hW/(kT T ) 
xu = No/Nn 

X n l = N H(nJ) /Nn 

x e = N e /Nn 
x = (y— u„)/(v n A K ) 
x = (v- v c )/(v c Au) 

X n = (Vn — Vc)/{v c A H ) 

y = h(v - ^ Lya )/(fcT r ) 
y D = h(w — WRjKkTr) 



A x = V2^ m /(m x e 2 ) 

Av e = I/LyqAe 



Eq. ([10 I 
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where Vhya ~ 2.47 x 10 Hz is the Lyman a frequency, 
and v c = I^Lya is the Lyman- limit frequency. The spon- 
taneous emission rate (Einstein A-coefhcient) in the Ly-n 
transition is IMl: 



2 13 7T 3 (l_ 

3 2 n 3 ~ 



2 13 7T 3 1 I/ 3 

n»i 3 2 exp(4) n 3 c 2 



l\ 2 ™-5 ,.3 

n I 7n 

1 \ 2n+5 r 2 ' 



(3) 



(4) 



where a is the fine structure constant and ao is the Bohr 
radius. 

We define the ratio 



rr, 
J CM 



3x ls 



(5) 



which is the equilibrium value of the photon occupation 
number at the Ly-n transition frequency. 

The Sobolev optical depth for the hydrogen Ly-n tran- 
sition is I2"2"ll551: 



3c 3 iy H xi s 



A 



3.T 



Is 



(6) 



where H(z) is the Hubble expansion rate. In all that 
follows, we will neglect stimulated emission in the Lyman 
lines, as the photon occupation number near Ly-n is of 
order f™ q <§; 1 (the largest value is for n — 2 and is less 
than 10 -11 for z < 1600). In particular, we can neglect 
the last tern J in the expression for the Sobolev optical 
depth Eq. 

Finally, we will refer to the matter temperature as T m 
and the radiation temperature as T r . In practice, the 
matter temperature is locked to the radiation tempera- 
ture through Thomson scattering, and the relative dif- 
ference between the two is below a percent until redshift 
z fa 500 El M- 



B. Line processes 

Consider an excited hydrogen atom in the np state; it 
has two mutually exclusive fates. 

The first possibility is that it reaches another excited 
state n's or n'd, with n' =/= 1, either through a spon- 
taneous or stimulated decay if n' < n or following the 
absorption of a CMB photon if n' > n. It can also be 
photoionized by a CMB photon. 

The second possibility is that the atom spontaneously 
decays to the ground state, emitting a Ly-n photon. In 
principle, this decay can also be stimulated. However, 



1 When dealing with the possibility of cosmological masers in 
Sec. |VI| we will of course explicitly account for this term and 
use the appropriate expression for the Sobolev optical depth in 
a general (not necessarily Lyman) transition. 



even when accounting for non thermal distortions to the 
radiation field, the photon occupation number at Lyman 
frequencies is extremely small, and stimulated emission 
in the Lyman lines can be neglected. 

The probabilities of these two complementary fates 
are denoted p™b>Psc respectively (the justification of the 
notation will become clearer in the next paragraph). 
Given the width (or inverse lifetime) T np of the np state 
(the sum of the rates of all transitions depopulating this 
state), they are given by: 



A 



1C 



np,ls 



np 



i 



Kb- 



(7) 



If we now assume that the considered atom was initially 
in the ground state and reached the np state after the 
absorption of a resonant Ly-n photon, the two fates men- 
tioned above can be described in a two-photon picture. 

First, if the atom reaches another excited state n' < n, 
the overall reaction H(ls) + 7(Ly-n) — > H(n'i) + 7' is 
a (possibly stimulated) Raman scattering event. If the 
atom absorbs a CMB photon and reaches a higher excited 
state n' > n (or gets photoionized), the overall reaction 
H(ls) + 7(Ly-n) + 7' -> R(n'l) (or -> e~ + p+) is a 
two-photon absorption (or two-photon photoionization) 
event. In these cases we will refer to the absorption of 
the Ly-n photon as a true absorption event, in the sense 
that the photon is destroyed in the process. The emission 
of a Ly-n photon following the inverse reaction chain will 
be referred to as a true emission event. 

Secondly, if the atom decays back to the ground state, 
the overall reaction H(ls) + 7(Ly-n) — > H(ls) + 7(Ly- 
n) is a Rayleigh scattering event (which in what follows 
we will refer to as a scattering event for short). In that 
case the incoming and outgoing photons have the same 
frequency in the atom's rest frame. Their frequencies in 
the comoving frame (frame in which the CMB appears 
isotropic) are Doppler-shifted with respect to the atom's 
rest frame frequencies. Since the Doppler shift depends 
on the relative orientation of the photon propagation di- 
rection and the atom's velocity, the frequencies of the 
incoming and outgoing photons in the comoving frame 
are in general different. They are however statistically 
correlated, as will be described in Section [II B 2| 



1. True absorption and emission 

The rate of true emission of resonant Ly-n photons at 
frequency v, per H atom, per frequency interval, is given 
by (2E1C3): 



\(n'#l),i 

x P n sc <p n {y). (8) 

In the above equation, R n ii^ np is the radiative transition 
rate per hydrogen atom from the n'l state to the np state, 
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including stimulated transitions, a np is the direct recom- 
bination coefficient to the np state, including stimulated 
recombinations, and ^Pniy) is the line profile, which has 
the Voigt shap^] 



where 



<fo(x;a) 



n^H 



+ 00 



tt 3 / 2 a 2 + (x- t) 

is the dimensionless Voigt profile, 

z T m 



(9) 



jdt, (10) 



A 



H 



IkTy. 



2.35 x 10" 



1 



m H c" \1100 T r 

is the dimensionless Doppler width and 



-1/2 



1 A 



np,ls 



4tw„A h pS c 47ri/ n A H 



(11) 



(12) 



is the dimensionless Voigt parameter of the line. 

For small Voigt parameters a <C 1, which is the case in 
all Lyman lines at the epoch of recombination, the Voigt 
profile has the well known asymptotic behaviors in the 
line center and in the damping wings: 



(f>v{x;a) 



4=e- a 



, x ^ x a 

, X X*, Xn 



where the transition scale x a is the solution of 



xle 



(14) 



In general, for a <C 1, x a ~ 3. 
Following Ref. [2S], we define 



(cm) 



3xi s Y n pp 



(15) 



a.b 



so that the true emission rate per H atom per frequency 
interval can be rewritten as: 



•^Icm =P"b 3x ls A n p )ls ^ n (l/)/ ( ™ m) . 



(16) 



The rate of true absorption of resonant photons is sim- 
ply the total rate of absorption times the true absorption 
probability. The absorption profile differs from the emis- 
sion profile by a factor e h ( v ~ Vn " ( kT ' , where T = T Y in the 
wings (because the low-energy photon of the two-photon 



process comes from a black-body distribution and the 
excited states of hydrogen are near Boltzmann equilib- 
rium with each other at temperature T r [J3J [S3]), and 
T = T m in the Doppler core, where atomic recoil tends 
to equilibrate the radiation field with the thermal veloc- 
ity distribution of the atoms. The rate of true absorption 
of resonant Ly-n photons at frequency v, per H atom, per 
frequency interval, is therefore given by: 



•^lab = -P2h^ x laA np ,u(p n (v)e hT " fw 



(17) 



The net (uncompensated) rate of true emission of reso- 
nant Ly-n photons at frequency v, per H atom, is there- 
fore: 



,ab 



p" b 3x ls A npAs ip n (v) 

/(em) ~~ e 



fu 



(18) 



It will be useful in what follows to cast this expression 
into a different but equivalent form, which can be done 
with the following considerations. 

The rate of change in the population of the np state is: 



Xnp ^ ^ x n 'iR n 'i—y n p -(- x e XpN^o. n p x n pV n p 



(em) 



xisR 



Is— >np 



^np^np- 



(19) 



The radiative rates (of order ~ 10 8 s -1 ) are many or- 
ders of magnitude larger than the overall recombination 
(13) rate, which is of order the Hubble rate H ~ 10~ 13 s _1 . 
The population of the excited states can therefore be ob- 
tained to an excellent accuracy by using the steady state 
approximation and setting x np = in the above equa- 



tion. Setting the left hand side of Eq. (191 to zero, we 
can solve for / ( ™ m) : 



rn 

J (ex 



1 



xi s R 



Is— ¥np\ 



(20) 



ab 



The total (including both true absorptions and absorp- 
tions leading to a scattering) Is — > np (forward only) 
excitation rate per H atom is given by: 



XlsRls 



¥np 



= 3x ls A npAs / ip n (v)e *r n f u dv. (21) 



This finally gives us the following relation for //" m j [13 



Pab/(em) 



fen-Pi / ¥>n(f)e kT " Udv. 



(22) 



We can now rewrite the net rate of true emission of reso- 
nant Ly-71 photons at frequency v, per H atom, in a form 
exactly equivalent to Eq. (18): 



2 The Voigt profile can be derived quantum-mechanically, in the 
two-photon picture, when one neglects the variation of multi- 
plicative factors v/u n across the line and uses the resonance ap- 
proximation. See e.g. Ref. [53) for fits to the correct profile in 
the case of the Lyman a line. 



+Psc 



, ab 



= 3xi s A nPt i s ip n (v)l /, 



eq 



e f v 



(23) 
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2. Coherent scattering 

The rate at which resonant scattering removes photons 
from the line, at frequency v, (in photons per second per 
H atom per frequency interval) is the total absorption 
rate times the scattering probability: 



K\ 



= -p" c 3x ls A npAs ip n {v)e hi "*T n) f u , (24) 



The rate at which resonant scattering injects photons in 
the line, at frequency v, depends on the absorption at all 
other frequencies, since incoming and outgoing photon 
frequencies are correlated: 

Pn {v\ V ') Vn { V ')e hi ^ 1 U,r\v' , (25) 

where p n (y\v') is the probability that the outgoing pho- 
ton has frequency v in the comoving frame given that 
the incoming photon had frequency v' . It accounts for 
the random thermal motions of the scattering atoms (at 
temperature T m ) and depends on the angular probability 
distribution of a scattering event. It it is normalized: 



for all v' ', / p n {v\v')dv = 1, 



and to respect detailed balance, must satisfy 

Pn{v\v')ip n (v') = p n (v'\v)ip n (v). 

The scattering redistribution kernel 



R n (v, v' 



--p n {v\v')ip n {v')& ( "^^ 



(26) 



(27) 



(28) 



is calculated in Ref. |57j (in which, however, atomic recoil 
during a scattering event is not accounted for). 

The most general form for the rate of change of the 
photon field through resonant scattering is given by 
Eqs. (24) and (25). However, in the case where the ra- 



diation field is smooth on the scale of a characteristic 



frequency shift in a scattering (A, 



,2x1/2 = 



2/ n A H , one 



can approximate the integral operator by a second-order 
differential operator. The rate of change of M v due to 
scattering is then given by a Fokker-Planck equation |51l - 
I53"II55H6"0] . accounting for scattering as a diffusive process 
in frequency space, with a systematic shift (or drift) due 
to recoil and Doppler boosting: 



FP 



' Psc3^1s^np,ls 

d_ 

dv 



kT„ 



(29) 



3. The radiative transfer equation 

In the vicinity of the Ly-n line, the photon occupation 
number evolves under the influence of the resonant pro- 
cesses described above, as well as eventual non resonant 



processes that may act in the vicinity of the line. The 
time-dependent radiative transfer equation (the Boltz- 
mann equation for the photon fluid) can be written in 
the general form: 



dt 



dt 



dfu 
dv 



,ab,; 



f»\ 



(30) 



where d/dt is the derivative along a photon trajectory, 
and /i/| nr groups all processes that are not resonant with 
the considered line (these could include absorption and 
emission from neighboring lines for example). The con- 
tribution of resonant processes is obtained with the con- 
version: 



f«\ 



„ab,i 



c^Vh 
8irv 2 



(31) 



Neglecting the variation of multiplicative factors vjv n 
across the line, and using the definition of r„, Eq. (JsJ) , 
as well as Eq. (23), the most general expression for the 



resonant term can be written 
1 



Hv 



f»\ 



,ab,i 



e kT 



h(u' - 

H> n {y)Wn(v )e " 



f - f r ' 
Rn(v,v')] U,dv'. (32) 



The radiative transfer equation, Eq. (30), with f u \ em>a ,b,sc 
given by Eq. (32), is therefore, in the general case, a 
time-dependent, partial integro-differential equation. As 
a result, it is computationally expensive to solve without 
further approximations. 

If the radiation field varies on a frequency scale large 
compared to a Doppler width, then we can use the 
Fokker-Planck operator for the scattering term, and the 
above term can be approximated by: 



■ — /J 

Tfy lem,ab,sc 



7v-/i 



(em) 



d_ 

'■ o 

ov 



v 2 A 



dv 



kT T , 



(33) 



More insight can be gained by considering some char- 
acteristic scales of the problem. Using the asymptotic 
expansion for the wings of the Voigt profile, Eq. ( 13 1, we 



obtain the total optical depth for true absorption in each 
damping wing: 



^ab,wing 



Pab T n- 



(34) 



If T at>,wmg ^ |- ne d am pi n g w ings are optically thin to 
true absorption, and one can use the Doppler core ap- 
proximation to the Voigt profile in the radiative transfer 
equation. If r„ ,wmg 3> 1, then the wings are optically 
thick for true absorption up to a detuning from the line 
center: 



Wr, 



; Pab T n fnA H > V n A H - 

7T 



(35) 
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In that case, the radiation field is near its equilibrium 
value f v « f(l m ^~ h ^~ Un ^' lkT within a detuning from the 
line center \v — v n \ < W n . 

One can similarly define the optical depth for resonant 
scattering in each damping wing, r^ c ' wing . If r^ c ' wing > 1, 
we can see using dimensional analysis that resonant scat- 
tering, when described with a Fokkcr-Planck operator, is 
effective up to a characteristic detuning from the line 
center 



S n = 



1/3 



(36) 



The physical meaning of this quantity can be understood 
as follows. A photon with initial frequency v can diffuse 
(in frequency space) to the line center in a characteristic 
time 



At 



diff 



(I/- 



1 



(^„A H ) 2 cx ls N-na n {vy 



(37) 



where a n {v) is the cross-section for resonant scattering: 

(38) 



3 c 2 

a n (v) =p^ c -—-^A np>ls ip n {v). 



The time it takes for the photon to redshift from v to v n 
(if v > u n , or from v n to v in the opposite case) is 



At 



redshift 



(39) 



Using the damping wings approximation for <p n {y) 
Eq. (131, the definition of the Sobolev optical depth 
Eq. (|6|, and Eq. (36), we obtain: 



Atdiff 



At 



redshift 



\v - v n \ 

Si, 



Therefore the radiation field will reach the equilibrium 
spectrum /„ oc e~' lI/ / fc T m within a detuning from the line 
center \u — u n \ < S n , due to the very fast redistribution 
of photon frequencies through resonant scattering. 

As an illustration, we show in Fig. [I] the parameters 
W2,<?2 for the Lyman-a line, extracted from the MLA 
code described in Ref. [44]. We see that for z > 800, 
S2 > W2 > ^LyaAu and at all relevant times S2 ^> 
i/LyaAjj. According to the above discussion, the radia- 
tion field in the vicinity of Lya is therefore smooth on a 
frequency scale Av ~ 6>2 around line center, and the use 
of the Fokker-Planck operator for resonant scattering is 
well justified. 



C. Net decay rate in the Lyman-n line 



<S 2 /iajA h . 
W2/1/2AH 



800 



1000 



1200 



1400 



1600 



FIG. 1. Characteristic widths W,S for the Lyman-a line, 
in units of Doppler widths, as a function of redshift for a 
standard recombination history. 



recombination, the quantities of interest are the net (un- 
compensated) bound-bound and bound-free rates, which 
are required to evolve the atomic level populations and 
the free electron fraction in a MLA code. In particular, 
the net rate of np — > Is decays is given by: 



•Enp— fls — %ls 



"np\ u 



dt 



J n,av = J a ; 1 



.ab 



(41) 



Using Eq. (23), we see that this can be written: 



J ec 



eq 



<p n (v)<lv, (42) 



(40) where we used the fact that the term proportional to p™ c 



in Eq. (23 1 integrates to zero. 



As can be seen in Eq. (42), the net decay rate in the 



Ly-n transition depends on the radiation field. The latter 
in turns depends on the atomic level populations through 
r n , p™ c and f" q , as well on eventual non-resonant processes 
acting in the vicinity of the line, as can be seen from the 
radiative transfer equation Eqs. (30), (32) or (33) . 



D. The Sobolev approximation 



Equations (30) and (32) give the most general forr 



of the radiative transfer equation, in a homogeneous ex- 
panding Universe. It has no analytic solution, because 
of its complexity and since it requires the values of the 
level populations, which in turn depend on the radiation 
field. It can however be simplified and decoupled under 



The exact shape of the radiation field in the vicinity of 
resonant lines is by itself of interest when predicting po- 
tentially observable spectral distortions to the black body 
spectrum of the CMB [37] . In the context of cosmological 



We neglected the variation of phase-space factors v/v n across the 
line, as well as stimulated emission and stimulated scatterings. 
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some approximations, for which an analytic solution ex- 
ists. The Sobolev approximation [22] [55] relies on the 
following assumptions: 

(i) No non-resonant proc esses act in the vicinity of the 
line (/„| nr = 0in Eq. g). 

(ii) Steady-state : the time it takes a photon to red- 
shift across the line is of order w/His n , where w is 
the characteristic line width. If the line is very nar- 
row, w -C v n , then this time is much smaller than 
the Hubble time, therefore physical quantities, such as 
nn,T T ,T m , H,x e ,x n i vary very little during the time it 
takes a photon to rcdshift through the line. As a conse- 
quence, one can neglect the time dependence in Eq. (30 1. 



(iii) Equal absorption and emission profile. This assump- 
tion also derives from the assumption of an infinitesimally 
thin line, in which case one can take the exponential fac- 
tors to be unity in Eq. (17). 



(iv) Complete redistribution of emitted photons. Mathe- 
matically, this means that p n (y\v') is independent of v 1 
in Eq. (25). Because of assumption (iii), it is also as- 
sumed that p n (p\v') = ip n (i>), i.e. the scattered photons 
are completely redistributed over the line profile. This 
implies that R n {v,u') = (p n (v)(p n (v') : and therefore the 
integral in Eq. (1321) vanishes (taking the exponential to 



be unity). 

The radiative transfer equation, under these assump- 
tions, becomes the much simpler ordinary differential 
equation: 



It has the analytic solution 

P+OO 

U = A + (/? - A) exp -T n / <p n (l/)dj/ 



(43) 



, (44) 



where /? is the photon occupation number at the blue 
side of the line. The photon occupation number at the 
red side of the line is therefore 



/- = /e"q + (/? - /e n c)e- 



Integrating Eq. (43) from — oo to +00 gives 



f n_ f n 



(45) 



(46) 



This, combined with Eqs. (|42j) and (45j finally gives us 
the standard Sobolev expression for t 
the line: 



re net decay rate in 

np^ls — 3xi s A„ P: i s Ps,n [fcq ~ /+] > (47) 



where 



Pi 



1-e" 



S,n 



(48) 



is the Sobolev escape probability. In the case of Lyman 
transitions, r n 3> 1, and Ps, n ~ l/r„. The net decay rate 
becomes, using Eq. {6}: 

= 8 ^ N f [A - fl] , (49) 



which is simply the rate at which distortion photons red- 
shift across the line. 

The Sobolev approximation provides relatively accu- 
rate net decay rates, despite the multiple assumptions 
that it relies on, and recombination histories currently 
used for CMB analysis [33] HH] > which use this approx- 
imation, are accurate at the percent level. The level of 
precision required by upcoming CMB experiments has 
motivated recent work to obtain more accurate solutions 
to the radiative transfer equation and net bound-bound 
rates in the optically thick Lyman lines. In this work we 
investigate previously ignored radiative transfer effects, 
and quantify as much as possible the errors made by the 
inevitable approximations that still need to be made. 



III. THE LYMAN ALPHA LINE 

The Lyman-a transition is one of the bottlenecks of 
hydrogen recombination. Electrons recombine to the ex- 
cited states of hydrogen, from which they eventually cas- 
cade down to the n = 2 state. They can then reach the 
ground state either by a two-photon decay from the 2s 
state, with rate A 2s i s rs 8.22 s _1 , or from the 2p state, 
by redshifting out of the Lyman-a resonance, with rate 
^bpis-Pesc, where P esc is the escape probability. Due to 
its substantial impact on the recombination history (see 
for example Fig. 11 in Ref. 07]), the net decay rate, or 
equivalently the escape probability in Lya has been stud- 
ied extensively, including time-dependent effects [47] [58] , 
two-photon processes gH [50], resonant scattering [51T - 
1551 I5TI 152"] , and Thomson scattering 52J . In this section 
we consider two additional effects: the non-local aspect 
of Thomson scattering (Section |lII A ), and a quantitative 
estimate of the effect of deuterium on hydrogen recombi- 
nation (Section MB). 



To simplify the notation we drop the subscripts and 
superscripts "2" in this section, and all the quantities 
previously defined implicitly refer to Lya. 



A. Thomson scattering in Lyman-a 

Thomson scattering in the vicinity of resonant lines 
was investigated in the context of helium recombination 
[313] with a Monte Carlo method and found to lead to 
negligible changes to the recombination history. Its effect 
on the hydrogen Lya line was investigated recently [52] 
and found to lead to negligible corrections to the escape 
probability and the recombination history. In Ref. [52 
however, electron scattering was described with the Kom- 
paneets equation, which is not valid in the context of cos- 
mological hydrogen recombination, as we argue below. 
Here we provide a more rigorous treatment of Thomson 
scattering, using the full redistribution kernel, which we 
incorporate in the Lyman-a radiative transfer code de- 
scribed in Ref. [53]. 
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The strength of Thomson scattering is characterized 
by its differential optical depth, flat in frequency in the 
non relativistic limit |3"01: 



N H x e a T c 

Hvhya 



(50) 



where ctt ~ 6.65 x 10 -25 cm 2 is the Thomson cross sec- 
tion. 

Thomson scatterings can affect the recombination his- 
tory if they take place within the characteristic width 
W over which the Lyman-a line is optically thick for 
true absorption. We show in Fig. [2] the mean number of 
Thomson scatterings within a detuning W of line cen- 
ter, rj e W '. We see that it peaks at ~ 0.08 for z ~ 1375, 
and remains above 0.001 for z > 1000, which suggests 
that Thomson scattering is potentially important at the 
subpercent level and should be carefully accounted for. 

The rate of change of the number of photons per unit 
frequency per hydrogen atom due to Thomson scattering, 
neglecting stimulated scatterings, is: 



-K 



(51) 

where Rt(v' — > v) is the electron scattering kernel. 

If the radiation field is smooth on the scale of a char- 
acteristic frequency shift during a scattering Ai/ C = 
v\jy a \/ IkTmj m e (? , then the integral operator for elec- 
tron scattering can be approximated by a Fokker-Planck 
operator, accounting for diffusion and drift in frequency 
space with rates [63] : 



d(Ai/ 2 ) 

dt 
d{Au) 
dt 



Nnx c a^c v 

NnX e <JTC v 



2 2kT m 
m e c 2 



hv 



rricC* 



(52) 
(53) 



The corresponding Fokker-Planck equation is known as 
the Kompaneets equation: 



\r l FP at kTm 
1 1 m c z 



d_ 
dv 



d 
dv 



h 

Wr 



(54) 



This is the approximation that was made in Ref. |52j . 
However, due to the small mass of the electron, the char- 
acteristic frequency shift during an electron scattering 
event Av e can be larger than the characteristic width 
over which the radiation field changes in the vicinity of 
the line. This characteristic width is of order S (which is 
> W at all times, see Fig. [T]), set by frequency diffusion 
due to resonant scattering near line center (see end of 
Section |II B 3[ ) . We see from Fig. [2] that Av e > S at all 
times, and therefore electron scattering cannot be consid- 
ered as a diffusive process and the Kompaneets equation 
is not valid in this context. 




FIG. 2. Top panel: characteristic number of electron scatter- 
ings within the region where the Lya line is optically thick 
for true absorption. Bottom panel: ratio of the characteris- 
tic frequency shift during a Thomson scattering event to the 
characteristic width over which the Lya line is smoothed out 
by frequent resonant scatterings. 



We have implemented the correct integral scattering 
kernel given by Eq. (51 1 in the Lyman-a transfer code 
developed in Ref. [S3] . Accurate approximate expressions 
of the electron scattering kernel Rr^(v' — > v) are given in 
Ref. [63] . For the purpose of our calculation, we only 



need the kernel calculated in the non-relativistic limit 



adequate here, for a dipolar angular distribution 
set 



R T {u' 



Rt{v v') = 



v — V 



We 



,(55) 



where the dimensionless kernel TZ is given bjj^] 63] 164] : 



1 



exp 



4 



3 + r 



1 

id 



P 4 



|/3 1 erfc 



2 



(56) 



We moreover require that detailed balance is satisfied, 
i.e. that the Planck spectrum M v cx v 2 Qr hv ^ kTT ^ (in the 
limit hv 3> kT m valid here) is preserved by imposing: 



Rt(u -> v 1 ) 
R T {v' v) 



exp 



(57) 



When evolving the number of photons per H nucleus per 
frequency bin in the i th frequency bin, iVj = M Vi ViA In v, 
the radiative transfer code uses a backward Euler method 
which requires inverting the matrix equation: 



Mij(t + At)Nj(t + At) = Ni(t). 



(58) 



4 There is a typo in Ref. Iti4l : erf should be erfc. 
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The matrix to be inverted, M, is tridiagonal in the case 
where only absorption, emission, and resonant scatter- 
ing (described by a Fokker-Planck operator) are present. 
Thomson scattering breaks this tridiagonality, which ren- 
ders the system prohibitively time-consuming to invert 
(M is a 801x801 matrix in our lowest resolution run). 
However, we can use the fact that Thomson scattering 
is only a perturbation to the radiative transfer equation. 
Therefore, M = Mo + <5M, where Mo is an easily invert- 
ible tridiagonal matrix, and the perturbation <5M due to 
Thomson scattering is such that its eigenvalues are al- 
ways small compared to those of Mo. We can therefore 
invert the perturbed matrix using the expansion: 
-l 



(M + SMy 1 = Mq - M i (,5M)M ( 7 i 
+ Mo X (5M)Mo 1 ((SM)M,7 : 



.(59) 



We find that the second order of the expansion is usually 
sufficient, with a maximum change of the net decay rate 
in the line of 1.5 x 10~ 5 between the first and second 
order. 

We show the resulting changes in the free electron frac- 
tion in Fig. [3] We can see that at early times, z > 1350, 
Thomson scattering delays recombination. Indeed, the 
relatively large frequency changes during electron scat- 
terings allow photons to be moved from the red side of 
the line to the blue side of the line, and vice versa. Be- 
cause of the large jump in photon occupation number 
across the line, the net photon flux is from the red side to 
the blue side. As a consequence, some escaping photons 
are reinjected into the line, where they can be absorbed, 
which decreases the escape rate and delays recombina- 
tion. At later times, this effect is not so important as 
the radiation profile becomes smoother (Af e ~ <S, see 
Fig. [2]). The systematic frequency loss during scatter- 
ing event due to electron recoil starts to dominate, and 
Thomson scattering helps photons escaping out of the 
line and speeds up recombination. 

For comparison, we have also implemented the Kom- 
paneets equation ( |54[ ), in a similar fashion as resonant 
scattering (see Ret". |53| for details on the implementa- 
tion). We can see that using the Kompaneets equation 
does not represent accurately the physics of Thomson 
scattering, as it cannot capture the large frequency shifts 
at early times. The error in the correction is of order 
the correction itself, and it has the wrong sign at early 
times. However, the basic conclusion reached in Ref. |52j 
remains valid: Thomson scattering can indeed be safely 
ignored during cosmic hydrogen recombination, since it 
leads to corrections to the ionization fraction of at most 
Ax e /x e - ±3 x 10" 5 . 



B. Interaction with the Deuterium Lya line 

1. Motivations 

A second radiative transfer effect associated with the 
Lyman-a transition is the interaction of the hydrogen 
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FIG. 3. Changes to the recombination history due to Thom- 
son scattering. 



and deuterium lines. Due to the slightly larger reduced 
mass of deuterium, the Lyman-a frequency in deuterium 
= ^Lya(D) is shifted to a higher frequency than that 
of hydrogen, v-r = i / Lya(H). The relative shift, to first 



order in m e /m p , and with mrj+ ~ 2m p , is: 



2m r 



2.7 x 10~ 4 . 



(60) 



This separation is ~ 10 times the Doppler width of 
the hydrogen line so the D Ly-a line center lies in 
the blue damping wing of the H Ly-a line. Despite 
the tiny fractional abundance of primordial deuterium 

x D = 2.87to;2i x 10_5 i 65 . ' the D L y-« line is stiU opti- 
cally thick during cosmological hydrogen recombination, 

T D,Lya ~ ^DTH.Lya ~ 10 2 — 10 4 . This has motivated 
the authors of Ref. [46 to consider the possible screen- 
ing of radiation incoming into the H Lya line by the 
optically thick, bluer D Lya line. However, due to the 
optically thick nature of the H Lya damping wings, the 
deuterium problem is more complicated than a simple 
Lya(D)— >Lya(H) feedback prescription, as pointed out 
in Ref. [IB]. In this section, we set up the problem of 
the interaction of the H and D Lya lines, and show that 
there is no significant effect on the recombination history. 



2. Spectral distortions caused by deuterium 

Our first step is to understand the physical mechanism 
of deuterium recombination. The rates of D and H re- 
combination are tied together via the charge-exchange 
reaction 



D+ +H(ls) O D(ls) +YP 



(61) 



which has a forward rate coefficient of order ~ 
10~ 9 cm 3 s -1 [SBJ; at recombination-era densities of ~ 
500 cm -3 this implies an equilibrium timescale of ~ 
2 x 10 6 s, i.e. six orders of magnitude shorter than the re- 
combination timescale itself. Thus to a very good approx- 
imation, the deuterium ionization fraction tracks that of 
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hydrogen: 



and write 



x D + 



x m 



XD(ls) ^H(ls) 



1 



(62) 



where the difference in ionization energies of deuterium 
and hydrogen, AE\/k rj 41 K, is small compared to the 
matter temperature during recombination. 

Like hydrogen, deuterium may recombine radiatively 
and reach the ground state either via 2s — > Is two-photon 
decay or Lya escape; the net rate of recombinations 
(H+D) is simply the sum of the two rates, and charge ex- 
change (Eq. [61 ) distributes the bound electrons between 
H and D. It is clear that the absolute rate of transitions to 
the ground state via the optically thin 2s — > Is transition 
will be much smaller for D than H because of its lower 
abundance; however for the Lya channel in principle the 
rates could be comparable because the Lya(D) optical 
depth is smaller than the Lya(H) optical depth by a fac- 
tor of ~ 2d , and hence the escape probability is enhanced 
by a factor of x^ 1 . However, since Lya(D) is located 
blueward of Lya(H), all photons emitted in Lya(D) will 
be re-absorbed in Lya(H); and all photons absorbed in 
the Lya(D) transition would have been absorbed anyway 
had the deuterium not been present. Thus the net (inte- 
grated over time) number of recombinations that proceed 
via Lya(D) escape is zero. However, the re-absorption of 
the D photons is not instantaneous, so at any given time 
the presence of deuterium causes an additional distor- 
tion of the radiation field (on top of the distortion to 
the blackbody spectrum due to H(2p) — > H(ls) decays 
in the blue wing of H Lya). The number of additional 
distortion photons per hydrogen atom is 



U 



&1TV 2 

c 3 N Vl 



(H+D) 



(63) 



(H) 

where fh is the photon occupation number when only 
hydrogen is present, and /^ H+D ^ is its value when the 
presence of deuterium is accounted for. 

Note that U(t) —> at both early times (because the 
radiation field is thermal independently of the presence 
of deuterium) and late times (since after recombination 
is over, there are no more Lya photons produced). 

The correction to the rate of formation of ground-state 
atoms through the (H+D) 2p — > Is channel is then equal 
to the rate of creation of distortion photons: 



s|Lya(D) 



-Ax 



e|Lya(D) 



u. 



(64) 



3. Analytic estimate for the number of spectral distortion 
photons 

We may now obtain an analytic estimate for U(t). Fol- 
lowing Hirata [U] and Hirata & Forbes [S3] , we make the 
change of variables 



v = m 



kT r 



-y, 



(65) 



fu 



( f H 



(66) 



where = x H(2p) / (3x H ^ 
sionless widths: 



We also define the dimcn- 



W : 
S- 
2/d : - 



hW 

kT T 
kT T 



(67) 
(68) 
(69) 



Taking T m ~ T r and using the damping wing approxi- 
mation for the Voigt profile, the steady-state radiative 
transfer equation in the blue wing of H Lya takes the 
following form in the absence of deuterium |53j : 



dy 



W 
!J 2 
_d_ 
dy 



d*(H) 



S 



dy 



*(H) 



(70) 



with boundary conditions ^( H )(+oo) = (neglecting 
Ly/3 — > Lya feedback in this problem for simplicity) and 
*( H )(0) = 1. The presence of deuterium modifies this 
equation by adding a term due to absorption and emis- 
sion by deuterium at y = y^: 



dy 



W 

y 2 

d_ 

dy 



(e^< H + D ) - 1 
d^( H +°) 



S 



^(H+D) 



+ T B w(y)(c y - VD ^ ii+D Hy)-yD) (71) 

where td = TD,Lya j ( ^d(?/) = h/kT T ipu{^) is the dimen- 
sionless emission profile in the deuterium line, and we 
have used 



*n = 



J cq 



k.Ty 



fH 

J eq 



(72) 



with /° = x D(2p )/ (3x D(ls) ) 



Several implicit assumptions went into Eq. ( 71 



First, we assumed complete redistribution in the deu- 
terium line (i.e. did not account for partial redistribution 
during resonant scattering events). This approximation 
is well justified since the deuterium line is dominated by 
its Doppler core: the damping wings are only marginally 
optically thick, and the differential optical depth in the 
deuterium wings is always much smaller than that of hy- 
drogen anyway (see for example Fig. 7 of Ref. [H]). In 
this case complete redistribution is a good approximation 
since both partial and complete redistribution have a sim- 
ilar characteristic frequency width, the Doppler width of 
the line. 
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In addition, we assumed that the Fokker-Planck rep- 
resentation of the scattering operator (for scatterings by 
hydrogen atoms) remained valid in the vicinity of the 
deuterium line. This can only be valid if the radiation 
field remains smooth on the scale of a Dopplcr width, 
even with the presence of the narrow optically thick deu- 
terium line. This is the case since (v^ — ^h) 3 <C S 3 at 
most times (see Fig. [I] with vd — i^h ~ 10 x vhAh). 
Therefore, according to the discussion in Section II B 3 



diffusion occurs on a much shorter timescale than red- 
shifting, and photons emitted in deuterium Lya decays 
at the deuterium resonant frequency are rapidly redis- 
tributed on both sides of the line by frequent scatterings 
by hydrogen atoms, smoothing out the radiation field in 
the vicinity of the D Lya line. 

Finally, one should in principle rescale the optical 
depth for hydrogen true absorption and scattering by a 
factor (1 — id), in both the blue and red wings (in par- 
ticular, change W and S to (1 — xu)W and (1 — xu)S 
in Eq. ( 71 )). However, Lya decays in the damping wings 
only lead to a ~ 1% correction to the recombination his- 
tory [44]. Therefore the correct rescaling would lead to 
a negligible 0(1O -7 ) correction, and we simply use W, S 
in Eq. ( 71 ) and neglect corrections to the decay rate in 



the red wing. 

Our last approximation will consist in taking tpr>(y) = 
3(y — yn)- Indeed, the radiation field will vary on a char- 
acteristic scale S, which is much larger than the width 
of the Doppler-broadened deuterium line (see Fig. [T] the 
Dopplcr width of the deuterium line is V2 times smaller 
than that of the hydrogen line) . The presence of the op- 
tically thick deuterium line in the blue wing of H Lya 
therefore amounts to imposing an additional boundary 
condition: 



*( H + D )(y D ) = *( H + D >( 2 ,+ ) = * I 



(73) 



where the continuity of ^( H + D ) across the deuterium line 
is insured by the very fast redistribution of emitted pho- 
tons in the frequency domain by the frequent hydrogen 
resonant scatterings. Eq. (71) is otherwise identical to 



Eq. ( 70 ) for y > y D and y < yn- We define 

A* D = * D -*( H )(?; D ). 
By linearity, we have 

^(H+D) _ ^(H) = x ^ 



(74) 



(75) 



where the function i/j(y) is the solution of the second 
order linear homogeneous equation 



dy 



W 

^e*V 



d 

dy 



S fdtp 



with boundary conditions 

^(+oo) = ^(0) = 



y 2 \dy 



(76) 



(77) 



The number of additional distortion photons per hydro- 
gen atom due the presence of deuterium therefore be- 
comes 



U 



c 3 N H h 



eq 



where 



AX(W,S, m ) 



+ OO 



ip(y)dy. 



(78) 



(79) 



Our last step is to estimate A^d ■ The population of the 
n = 2 shell of deuterium is controlled by recombinations 
to the excited states and Lyman-a decays (as mentioned 
earlier, decays to the ground state through the Lyman-a 
channel are dominant over two-photon decays from the 
2s state because of the relatively low optical depth of 
the D Lya line). We can obtain the contribution of the 
former with a Peebles-like estimate [2~4l: 



XB(2p) I 



-a B (T m ) 

1 / m c kT n 
3 V 2irh 2 



Nhx c x d + 

3 / 2 _„„ 



x ~D(2p) 



(80) 



where aB(T m ) is the case-B recombination coefficient, 
for which we use the fitting function given in Ref. |67j . 
The contribution from the Lya decays, Eq. (42), can be 
rewritten after some manipulations as: 

x J r D (^ D (y) (ey-y»¥ n+I >)(y) - tf D ) dy. (81) 



Using Eq. (71 1, and taking the limit <pv(y) — > 5(y — i/d) 
(i.e. assuming tfr>(jj) has support only in y^, ± e, and 
taking the limit e — > 0) we obtain: 



(y) (ey-y»¥ H+D Hy)-y D ) = 



s 



d*( H + D ) 



dy 



d$( H + D ) 



yi 



dy 



(82) 



where we used the continuity of vI>( H+D ) across the line 
(so the first derivative integrates to vE»( H+D ) (j/ D + e) — 
\jf( H + D )(y D _ e ) 0), and its finiteness (so the inte- 
gral of the term proportional to W vanishes as 6 ^ 0). 
The derivative of *( H+D ) may however be discontinu- 
ous across the deuterium line because of the presence 
of the delta function. Eq. (82) simply states that the 



net rate of decays in the deuterium line is balanced, in 
steady-state, by the flux of photons due to frequency 
diffusion across v^. Since there is no discontinuity of 
the derivative of \I/( H ) at j/d, we can replace v[r( H + D ) by 
$(H+D) _ ^(H) = A ^ D x V> in Eq. Our final expres- 



sion for the net rate of Lya absorptions is therefore 



£D(2p) 



8rHvi 



Lya 



c 3 iY, 



H V° q 



A* D AS, (83) 
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where 



AE(W,S, m ) 



s 


'dip 




dtp 






dy 




dy 





> 0. 



(84) 



Transition rates into and out of the 2p state of deuterium 
are many orders of magnitude larger than the overall re- 
combination rate which is of order of the Hubble rate. 
The population of the 2p state, and therefore A 'I'd, can 
be solved for in the steady-state approximation, setting 



« ±D(2p) = ±D(2p) Lcc + *D(2p) | Lyc 



(85) 



We first express £d(2p) in terms of A 'I'd, using Eqs. (66 1, 
(Tfl and ([741): 



x D(2p) 



= 3a; 



D(ls) 



*W(y D ) + A*! 



(86) 



We also need to compute \I/' h '(j/d), AH and AT numer- 
ically as a function of W, S, y-Q . We do so by solving 
Eqs. (70) and ( [76]) with a shooting method (see Ap- 
pendix A of Ret. |53|). Using Eq. (62) for ir D + and re- 
placing Xrj( 2 p) by the above expression in Eq. (80 1, we 
can then solve for A'I'd for a given "order zero re- 
combination history (i.e. computed without the deu- 
terium correction) and the corresponding parameters 
W(z),S(z), f^ q (z). We use the output of our EMLA code 
[IT] to generate such a history. 



4- Results and discussion 

We obtain a very small relative distorsion A'I'd, with 
a maximum value of ~ 10~ 4 at z ~ 1600, falling below 
10 -6 for z < 1400. As a result, the number of distortion 



photons per hydrogen atom, U, given by Eq. ( 78 ) is also 
extremely small, peaking at ~ 3 x 10 -9 at z ~ 1500, and 
D Lya decays have virtually no effect on the recombina- 
tion dynamics. 

To better understand the magnitude of the effect, it is 
more enlightening at this point to write some approxi- 
mate expressions in the case relevant here. 

We start by rewriting the net decay rate in the deu- 
terium Lya line, Eq. (83), in a more familiar form, u sing 
the definitions for A* D , V D and * H (y) Eqs. ([74), (f72| 



and (66 1 



;C D(2p)| Lyo 



-AS 



3c 3 A^ H a;D( 



X (xo(2p) - 3^(15)/^) 



(87) 



where /+> = Without the factor AH, this is just 

the usual Sobolev decay rate, Eq. (49) which is the rate 



at which photons emitted in the line redshift across the 
line. 

Now, as can be seen from Fig.[T] at most relevant times, 
W 3 <C S. Moreover, as the separation of the deuterium 



and hydrogen lines is approximately 10 Doppler widths, 
we also have y D <C S. In these limits, one can show that 
AH takes on the approximate form: 



AH 



Vd 



S 



» l. 



Using Eq. (40), we see that As ~ Ai re dshift/Atdiff- The 



meaning of Eq. (87) is now clearer: the net D 2p — > Is 
decay rate is given by the rate at which photons diffuse 
out of the deuterium line due to scattering by hydrogen 
atoms, rather than the rate at which they redshift out 
of the line as would be prescribed by the Sobolev ap- 
proximation. The fast frequency diffusion rate near vy> 
therefore enhances the decay rate in the line. Note that 
the decay rate per deuterium atom is already enhanced 
by a factor of i^ 1 relative to the equivalent rate in the 
hydrogen line, per hydrogen atom. Therefore the very 
fast 2p Is decays bring the deuterium 2p to Is ratio 
towards near equilibrium with the local radiation field, 
(if it were not for the distortions due to the optical 
thickness of the blue wing of hydrogen, the 2p state would 
therefore rather be near Boltzmann equilibrium with the 
ground state than near Saha equilibrium as it is the case 
for hydrogen). 

In the net recombination rate to the D(2p) state, 
Eq. (80), we can therefore approximate :ed(2p) ~ 
3s D (ii)7E, ~ e~ h ^~^/ kT 'x u x H(2p) . The second ap- 
proximation stems from the fact that (vy> — vn) 3 <C S 3 . 



With this approximation, and using Eq. (62 1, we obtain: 

(89) 



X D(2p)| 



£D^H(2p)| 



We also rewrite the net decay rate in the D Lya line, 
Eq. (83), as follows: 



x D(2p) 



Lya 



3; H(2p)| LyQ 



x A* D AH 



-i H ( 2 p)| rec x A4- D As, 



(90) 



where in the second line we used the steady-state approx- 
imation for the population of the H(2p) state. Setting 
a: D(2p) ~ 0, we therefore obtain: 



A*r 



Xu 

AH' 



(91) 



Now, roughly speaking, for j/d <C W, S, the integral 
AX = J ip(y)dy is of the same order of magnitude as 
the integral J 'i> li (y)dy (this can be checked numerically). 
This is only very approximate but allows us to rewrite 



c 3 N H h 



( f H 



AT • 



(92) 



where x+ is the number of distortion photons per hydro- 
gen atoms in the blue wing of Lya, due to two-photon 
transitions and diffusion in hydrogen (see Eq. (69) of 
Ref. [S3]). 
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We can now finally obtain an approximate expression 
for the number of distortion photons due to the presence 
of deuterium: 



U 



(93) 



The distortion in the blue wing of Lyman-a due to hy- 
drogen two-photon decays and diffusion itself leads to 
0(1%) corrections to the recombination history |4"4l 153] . 
The distortion due to the presence of deuterium is much 
smaller due to the very fast D Lya decays, which bring 
the D(2p)/D(ls) ratio very close to equilibrium with the 
radiation field at the frequency isq. As explained above, 
this fast decay rate is due to a relatively small optical 
depth in the deuterium Lyman line (smaller by a factor 
2d than that of hydrogen), and is further enhanced by 
the large diffusion rate near z/d, translating to the large 
AS in the denominator of Eq. (93). 

We therefore conclude that the effect of deuterium on 
the recombination history is at most 0(1CP 7 ) (due to the 
rescaling of the optical depth of hydrogen absorption and 
scattering in the damping wings, changing x + by factors 
of order id). Note that accounting for the presence of 
deuterium also induces changes to the background expan- 
sion of the order of xb ~ 10~ 5 . Moreover, depending on 
how the ionization fraction is normalized (with respect to 
total H+D of just H), there can be ambiguities of ~ 10 -5 
in x e . These modifications to the recombination histories 
can however be safely ignored given the level of accuracy 
required by Planck. 



IV. HIGHER ORDER, NON-OVERLAPPING 
LYMAN LINES (2 < n < 23) 

A. List of efficient Lyman transitions 

One of the strengths of the EMLA formulation |5T] is 
that it only requires the computations of effective tran- 
sition rates for a small set of "interface" states. These 
"interface" states should in principle be 2s, 2p, and all the 
higher lying p states. However, as exposed in Ref. [4"T] . 
only the lowest few Lyman transitions significantly con- 
tribute to the overall recombination rate, and one can 
neglect higher order Lyman transitions without loss of 
accuracy. The aim of this section is to verify this state- 
ment quantitatively. 

For that purpose, we computed a series of recombina- 
tion histories with RecSparse [3S], with n max = 30, arti- 
ficially setting the Lyman transition rates to zero above 
a fiducial transition, and compared them to the exact 
calculation where all Lyman transitions are included up 
to n max . Both calculations included feedback between 
neighboring Lyman lines with non-vanishing transition 
rates, and in both cases the photon occupation number 
incoming on the highest line considered was taken to be a 
blackbody. We find that neglecting Ly7 and above leads 
to changes to the recombination history \Ax e /x e \ of at 



most 8 x 10~ 5 , neglecting Ly5 and above leads to changes 
of at most 10~ 5 , and neglecting Lye and above leads to 
changes of at most 3 x 10~ 6 . We therefore conclude that 
it is sufficient to only include the Lya, /3, and 7 transi- 
tions when computing recombination histories accurate 
to the level of 0.01%. 



B. Resonant scattering in the low lying Lyman 
lines 



The Lyman-a line has a particularly high resonant 
scattering probability, as atoms in the 2p state can 
leave the state only through spontaneous Ly-a decay, 
or through absorption of a CMB photon. There is no 
lower energy state to spontaneously decay to other than 
the ground state. For higher order Lyman transitions, 
though, there are multiple allowed channels out of the 
np state, np — > n'l, with 1 < n' < n and I = 0, 2. There- 
fore the scattering probability defined in Eq. |7]) becomes 
comparable to the absorption probability, p" c ~ p" h (in 
Ly/3 * " ' ' 



'ail 



0.12 and p^ 1 



J ab 



0.16). From Eqs. (351, 



(361, we can see that in that case 



(94) 



The wings of the Ly/3 and Ly7 lines are optically thick 
for true absorption (this is the case for all Ly-n lines with 
n < 13). This implies (see Eq. (351 and corresponding 
discussion): 



1 < 



(95) 



where the the second inequality is a consequence of the 
first one. From Eq. ( 94 ) , we therefore obtain that for the 



low-lying Lyman lines above Lya, (S n /W n ) <C 1- This 
means that frequency diffusion is efficient only on a small 
fraction of the width over which the line is optically thick 
to true absorption. 

There is not a simple relationship between the ratio 
(S n /W n ) 3 and the impact of frequency diffusion on the 
Ly-n decay rate. However, the effect of frequency diffu- 
sion clearly increases with this ratio. Radiative transfer 
computations including frequency diffusion were carried 
out for the Lyman-a line, which showed that frequency 
diffusion leads to corrections of a few percents to the 
net decay rate in Lya [5H 153]. As shown in Fig. [I] 
(S2/W2) 3 3> 1 at most times for the Lya line, in con- 
trast with what we just showed for the higher-order lines. 
Therefore, we can expect that frequency diffusion would 
lead to corrections of much less than a percent to the net 
decay rate in the Ly/3 and Ly7 lines. Since Ly/3 decays 
themselves contribute of the order of a percent only to 
the overall recombination rate, sub-percent corrections 
to their rate can therefore be safely neglected at the level 
of accuracy required. 
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As a conclusion, resonant scattering in Ly/3 and higher 
order Lyman lines does not affect the recombination dy- 
namics to any significant level. 



V. OVERLAP OF THE HIGH-LYING LYMAN 
LINES (n > 24) 

A. Motivations 

In the previous section we have shown that Lyman 
transitions above Ly7 do not affect the recombination 
history to a significant level and can be ignored. How- 
ever, this relied on using the Sobolev approximation for 
the net decay rate out of the optically thick Lyman lines. 
As we show below, high-lying Lyman lines overlap with 
each other, or even with the continuum, which breaks as- 
sumption (i) of the Sobolev approximation that each line 
is isolated. In this section we address the consequences 
of this feature. 

For n > 13, the Lyman lines are dominated by their 
Doppler core (the damping wings become optically thin) . 
The condition for the two neighboring Ly-(n 4- 1) and 
Ly-n lines to overlap is therefore that their separation 
v n +i — v n ~ 2v c /n 3 becomes of order the Doppler width 
of an individual line, t'nAfj ~ ^^h- Overlap therefore 
occurs for 



> 



2 

Ah 



1/3 



44 



( T m 
\ 3000 K 



-1/6 



(96) 



For even higher order lines, the separation with the con- 
tinuum v c —v n = v c /n 2 becomes of the order of a Doppler 
width. This occurs for 



> 



= A~ 1/2 



206 



3000K ) 



-1/4 



(97) 



Intuitively, one can expect that line overlap amounts 
to adding new transitions R^p\ n i p between high-lying p 
states, as photons emitted in the np — > Is transition can 
be re-absorbed immediately in a neighboring Is — > n'p 
transition. More importantly, overlap with the contin- 
uum provides an additional recombination pathway. Di- 
rect recombinations to the ground state are usually con- 
sidered as highly inefficient as the resulting emitted pho- 
tons can immediately ionize neutral hydrogen atoms in 
their ground state. Chluba & Sunyaev [35] have con- 
sidered the possibility of continuum escape (similar to 
Lyman-a escape) , and have shown that it leads to negligi- 
ble corrections to the ionization history Ax e /x e ~ 10~ 6 . 
If overlap of the highest Lyman lines with the continuum 
is accounted for, it becomes possible for free electrons and 
protons to successfully recombine to the ground state of 
hydrogen, if the emitted photon subsequently excites an- 
other atom to a high- lying p state rather that ionizing it. 
The recombination event e~ + p — > Is + 7 immediately 
followed by the absorption event Is + 7 — > np thus cor- 
responds to an additional, indirect recombination event 



to the np state, to which a coefficient a4p can be asso- 
ciated. The reverse process corresponds to an additional 
photoionization rate, • 

Modern recombination codes account for the excited 
states of hydrogen up to extremely high principal quan- 
tum number n max > 200 [3H1 HO]: and recently up to 



500 



It is therefore important to quantify 



the impact of line overlap on the recombination history. 
In what follows we develop a formalism that generalizes 
the Sobolev escape probability method, and accounts for 
the overlap of the high lying Lyman lines. 

We start by evaluating the effect of photoionization 
and recombinations from and to the ground state on the 
radiation field and providing the relevant equations. 



B. Photoionization and recombination from and to 
the ground state 

The frequency-dependent photoionization cross- 
section from the ground state is, in the atom's rest 
frame: 



o"0) = (T g 



where 



9^2 



2 9 ^ 



Co 



-aa. 



I k 6.3 x 10" 18 cm 2 



(98) 



(99) 



3exp(4) 

is the photoionization cross-section at threshold and the 
function g(n) is such that g(n < 1) = 0, g(l) — 1, and 
varies on a scale Ak ~ 1 for k > 1 [54] . The Doppler- 
averaged cross-section is therefore: 



a{v;T u 



<t(i/[1 — mAh]) — -^du 



= cr 



9 



-[1-uAh] 



du(100) 



In what follows we use the distance to the Lyman limit 
in Doppler width units: 

(101) 



v c A u 



The Doppler-averaged cross-section can be rewritten as 
a{v;T m ) = 



1+XAH 5 ([1 + xAh][1-uA h ])- 



(102) 



For frequencies within a few Doppler widths from the Ly- 
man limit (| a; I ~ a few), the argument of g in the integral 
is 1 + O(Ah). Since g varies very little on a scale Ah, 
we can set it to its threshold value in the integral, g ~ 1. 
We therefore obtain the Doppler-averaged photoioniza- 
tion cross section, for frequencies within a few Doppler 
widths of the ionization threshold: 



cr(i/;T m ) RJ er 



1 - ^erfc(x) 



a Mx), (103) 
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where the last equality defines the dimensionless profile 

<t> c (x). 

The rate of photoionizations from the ground state per 
hydrogen atom per frequency interval (for the photoion- 
izing photon) is then: 



M 



phot 



-ca(v;T m )NH_x ls J\f„. 



(104) 



Detailed balance considerations show that the differen- 
tial recombination coefficient to the ground state, per 
frequency interval for the outgoing photon (in units of 
cm s _1 Hz _1 ) is given by: 



dais 
dv 



8irv 2 _ _Mi 
— ir—a(v\ T m )e k 



27rft 2 
m P kT„ 



.3/2 



(105) 



The rate of recombinations to the ground state per hydro- 
gen atom per frequency interval (for the emitted photon) 
is then 



M\ =N K . 



^dau 



(106) 



We define the continuum equilibrium occupation number: 
- 2 ' 2irh 2 ^ 3/2 



eq 



N H x j 



(107) 



The rate of change of the photon occupation number due 
to photoionizations and recombinations from and to the 
ground state, neglecting stimulated recombinations, can 
then be written as: 



f v \= -ca(v; T m )N H x ls f v -e / c c , 



.(108) 



C. The radiative transfer equation in the presence 
of multiple overlapping lines, and photoionization 
and recombination from and to the ground state 

The time-dependent radiative transfer equation for the 



photon occupation number /„, Eqs. (30), (32 1, in the 



presence of multiple lines, is (if one line is considered 
as the fiducial line, then the other lines and continuum 
absorption and emission constitute the "non-resonant" 
processes in that equation): 



Hp dt 



dv 



e fcr n , 



- r _ rn 
Ju .Icq 



Psc T n 



Hv 



fu 



- re 

J (_)• 



fvdv' 



(109) 



Since line overlap is expected to be a small correction to 
the recombination history, we can neglect "corrections to 
the correction" and make some approximations to sim- 
plify the calculations. As the potentially important effect 



is the violation of assumption (i) of the Sobolev approx- 
imation, we will lift this assumption, but keep the other 
three assumptions on which it relies, as we justify below. 

First, as we are considering the radiation field over 
a characteristic frequency width a few Doppler widths, 
which corresponds to changes in the scale factor Aa/a ~ 
10~ 4 , we can make the usual steady-state approximation 
by neglecting the time derivative. We also approximate 
the exponentials by unity since their exponents are of or- 
der ~ Aft-p^r- ~ 10~ 3 . In addition, we approximate the 
Doppler width of the Ly-n line v n An w ^ c Ah, neglecting 
corrections of order 0(n~ 2 ). Finally, we will assume com- 
plete redistribution for resonant scattering. The validity 
of the latter approximation is more difficult to precisely 
quantify, but it can be justified with the following ar- 
guments. Firstly, the resonant scattering probability p™ 
rapidly decreases as n increases due to the abundance 
of low-energy photons that can easily photoionize atoms 
in the np state or cause transitions to neighboring ex- 
cited states. We find, for T m = 3000 K, = 0.44, 
j4° = 0.28, ?4 00 = 0.17 and = 0.09. Secondly, for a 
Doppler-dominated line, partial redistribution is close to 
complete redistribution, in the sense that both distribu- 
tions have a similar characteristic width, of the order of a 
Doppler width. This contrasts with the Lya line, where 
complete redistribution can change photon frequencies 
by many Doppler widths due to the optical thickness of 
the Damping wings, and the distinction between the two 
types of redistribution is important. 

We work with the dimensionless frequency x defined 
in Eq. (101 1. We define the optical depth for continuum 



absorption, per unit Doppler width: 

§c{x) = T c (j> c (x), 

ca N H x ls 



where 



H 



We further define: 



1 



frAt 



$„(x) = T n cf> n (x) = r n —=cr 
$(jc) = ^ + $ c (a;) 



(110) 
(111) 

(112) 
(113) 
(114) 



The steady-state radiative transfer equation for f(x) be- 
comes, with the approximations justified above: 

~=J2 ~ /eq] +*c(s) [f(x) ~ / c C q ] (H5) 

n 

Since the continuum is optically thick, we set the bound- 
ary condition to /(+oo) = f^ q . Note that our treatment 
does not allow for any continuum escape in the absence 
of high-lying Lyman lines (the total optical depth for 
continuum absorption is infinite in our approximation). 
However, this has been shown to be negligible and lead 
to corrections to the ionization history Ax e /x e ~ 10~ 6 
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D. Generalized escape probability formalism 

1. Preliminaries 

Let us consider the probability distribution II„ of pho- 
tons injected with a total rate ri n j and a profile pinfa) 
(normalized to unity), that then undergo line and con- 
tinuum absorption. The evolution of the probability dis- 
tribution is given by an equation similar to the radiative 
transfer equation Eq. ( 109 1 (assuming complete redistri- 



bution and taking the exponential terms to unity): 



dn 



Hv 



At dt dv 



"inj^injM 

- ca{v)N n x ls Ti v . (116) 



The rate at which these photons are absorbed in the 
Ly-n transition is then Hvr n J </J n (^)II„d^, and the 
rate at which they are absorbed by the continuum is 
/ ccj(v)NyiXi s TI v <1v . Finally, these photons may also es- 
cape the set of overlapping lines by redshifting below 
their resonant frequencies. If we only consider lines above 
a given frequency v\ ow (below which lines can be con- 
sidered as isolated), then the escape rate is Hv\ ov ,H v ^ . 
In practice, v\ ow is many Doppler widths below v c (so 
^low *C —1), but still close enough to v c that we can 
approximate the escape rate by Hv c W v (x — ¥ -co). 

Now in steady state, the sum of all these rates must 
equal the injection rate (there are no other possible fates 
than those described above for the injected phot ons) . 
This can be checked explicitly by integrating Eq. (116) 
from v\ ow to +oo in the steady-state limit (with bound- 
ary condition n(+oo) = since no photons are injected 
at infinity). Therefore, the steady-state probability that 
injected photons are absorbed in the Ly-n line is 

P(inj n) = ~r„ f Vn (y)U v dv, (117) 

1 inj J 

the probability that they cause a photoionization is 



P(inj 



1 



ca(v)NnXi s Il v dv, 



(118) 



and the probability that the escape without being ab- 
sorbed is 



P(inj 



Hv, 



II„(a; — > -co). 



(119) 



In what follows, we apply this idea to develop a formalism 
for interline transition probabilities. 



2. Interline transition probabilities 



Taking advantage of the linearity of Eq. (115), we can 
decompose f(x) on a set of basis functions: 

/+oo 
My)G(x,y)dy, (120) 
-oo 



where the functions Vj(x) satisfy the linear inhomoge- 
neous differential equations: 



dx 



$>(x)vj{x) - (j>j(x), 



(121) 



with boundary conditions Vj{+od) = 0, and the Green's 
function G(x, y) satisfy a similar equation: 

dG 



d.r 



= $(x)G(x,y)-5(x-y), 



(122) 



with boundary condition G(+oo, y) = 0. The asymptotic 
behavior of G(x,y) at x 3> 1 is: 



G(x,y) 







,T a {x-y) 



x > y, 

x < y, 



(123) 



where we used $(+oo) = $ c (+oo) = t c . We therefore 
recover the appropriate boundary condition for f(x) at 
+oo: 



f(x » 1) w 1% 



+oo 



dy = 11 



(124) 



We can now use the results from the previous section. 
Using n„ = Vj(x)/(v c An), <p- m j M = <t>j{x)/{u c An), and 



H/Ah, we see that Eq. (121) is the steady-state 



version of Eq. ( |116[ ). Therefore the steady-state proba- 
bilities that a photon emitted in the Ly-j line is later ab- 
sorbed in a Ly-i transition, or subsequently photo-ionizes 
an atom in its ground state are, respectively: 



3>i(x)vj(x)dx, 



P(j -> c) = / $ c (x)vj(x)dx. 



(125) 
(126) 



The probability that a photon emitted in the Ly-j line 
escapes at x = — oo without being reabsorbed in any line 
or causing a photoionization is given by 



z) = Vj(-Oo). 



(127) 



Clearly, the probability of escape from the whole set of 
overlapping lines is vanishingly small, except possibly 
for photons emitted from the lowest lying line consid- 
ered as "overlapping" . The region of line overlap blends 
smoothly into the region of quasi-instantaneous feedback 
between neighboring Lyman lines, for n < 20—30. There- 
fore, even for photons emitted from the lowest line con- 
sidered as "overlapping" with the next higher line, there 
is still a near-unity probability of being reabsorbed quasi 
instantaneously in the next lower transition. Therefore, 
in practice, we have P(j — > esc) = for all lines consid- 
ered as overlapping, or simply close enough that feedback 
is quasi-instantaneous. 

Integrating Eq. (121) from — oo to +oo, we can see that 



these probabilities are complementary, as they should: 



^P(j^i)+P(j 



-> c) = 1. 



(128) 
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The Green's function G(x, y) can similarly be inter- 
preted as the steady-state number distribution for con- 
tinuum photons initially injected at the frequency y. The 
probabilities that a photon emitted at frequency y is ab- 
sorbed in the Ly-i transition, absorbed by the continuum, 
or escapes at -co are, respectively: 



P(y^c) = J ^ c (x)G(x,y)da 
P(y — > esc) = G(—oo,y). 



(129) 

(130) 
(131) 



Again, the probability for photons emitted in the con- 
tinuum to escape the whole set of high-lying overlapping 
lines is vanishingly small, which means that in practice 
we have P(y — > esc) = 0. 

Finally, one can check that these probabilities are in- 
deed complementary by integrating Eq. ( 122 ) between 
— oo and +oo: 



P(y -> c) = 1. 



(132) 



3. Net decay rate in the Ly-i transition 



From Eq. (42), the net decay rate in the Ly— i line is 



x ip ^ ls = 3x ls A ip , ls I [fl q - f(x)] (pi(x)dx 



f(x)$i(x)dx. (133) 



Using the expansion (1201 for f(x) and the definitions 
for the interline transition probabilities, we rewrite 



E/ TjXj p P(j 



A- i /"+ 00 



3 



We notice that Ai p ^ s /ri x tj = Aj P As ( approximating 
Vi w Vj w v c ). Also, using Eqs. (105) and (111), we 
obtain that 



3x ls Ai ^f^ c (y) = Nnx 2 e d ^ 



where = ^Ah^j 11 is the differential recombination 

ax *- 11 di> 

coefficient to the ground state per unit Doppler width. 
Using the complementarity relation (128), we can now 



rewrite the net decay rate in the Ly-i transition as: 



x 



ip— >ls 



X, 



ip 



?( ov ) _|_ fl( ov ) 
ip—t-jp ^ip 



E 



E4 

x ■ B [ov) - Nv,x 2 a {ov} 

^JP jp— up Iy tl^ e { -*-tp J 



(136) 



where we have defined the overlap-induced transition 
rates 



R 



(ov) 



Ai PiU P(i ->• j), 



(137) 



and the overlap-induced recombination coefficients and 
photoionization rates: 



Q 



(ov) _ f da-i 



(138) 



/3g v) = A iPtls P(i -> c). (139) 



4- Net rate of recombinations to the ground state 

The net rate of recombinations to the ground state, 
per hydrogen atom, is 

8nv? 



x u I a(u) 



■ h(v-v c ) 



dv 



x ls v c A H o- a / <j> c {x) [/ c c q - f(x)] dx, (140) 



where in the second line we took e feT m w 1 in the 
vicinity of the Lyman limit. Using again the decomposi- 
tion ( 120 ) and the definitions of the interline transition 



probabilities, we get: 



Xc— lis 



57TI/ 

3c 2 

where we have defined 



x ls v c A H a fe q x T c 

E r «- ip u ( 141 ) 



r. 



<j) c {x) 



1 



<t> c (y)G(x,y)dy 



dx 



dx 4> c (x) / dy 



8{y-x) - / & c (y)G{x,y) 



dy J dx <f> c (x) 

4> c (y) [i - P(y -> c)] dy 



5(1/ -x) - / $ c (y)G(x,y) 



(142) 



(135) 

where we have used the definition (|130[) of P(y — > c) in 



the fourth line, and in the last line we have used the 



complementarity relation (132). 



After some algebraic manipulations, we can cast the 
net rate of recombinations to the ground state in the 
following form: 



Xc— >ls — 



2 (ov) _ o(ov) 
X e Otjp XjpPjp 



(143) 



where the overlap- induced recombination coefficients and 
photoionization rates have been defined in Eqs. ( 138 ) and 



(139) 
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d t / d x (T Ly a = 6E8, T= 3000 K) 




-60 -40 -20 



FIG. 4. Differential optical depth per unit Doppler width, 
$0), for fiducial values r LyQ = 6 x 10 8 and T m = 3000 K. 
The dashed line is the contribution from Lyman lines, the 
dot-dashed line is the contribution from photoionization from 
the ground state, and the solid line is their sum. The position 
of the Lyman-30 to 100 resonant frequencies (spaced by An = 
10) is shown. 



5. Rate of change of the ground state population 



Adding Eqs. (136) and (143), we see that the net rate 
of change for the ground state population due to de- 
cays from the overlapping lines and recombination to the 
ground state is zero: 



X 



jp-*ls 



+ X c 



0. 



(144) 



0, ov 



This can be understood intuitively since any photon 
emitted from one of the overlapping lines is bound to 
be reabsorbed in a neighboring line or photoionize an 
atom from the ground state; similarly, a photon emit- 
ted after a recombination to the ground state will almost 
certainly be reabsorbed in a high-lying Lyman line, or 
cause a subsequent photoionization. In other words, no 
escape is possible from the series of high-lying lines, and 
the net decay rate to the ground state therefore vanishes. 
Line overlap can therefore only influence the recombina- 
tion history indirectly, through changing the populations 
of the excited states. 

We now turn to the numerical evaluation of the 
overlap-induced transition rates. 



E. Evaluation of the overlap-induced transition 
rates 

We show in Fig. [4] the optical depth per unit Doppler 
width, $(x). It can be seen that for x > x riov , is 
nearly constant. This can be explicitly derived as follows. 



From Eqs. 1^, pi), ([99]) and ( |111[ ), we can show that 

2 

T « ~ n a 3 r c- ( 145 ) 

n>l A H Ii' i 

Moreover, the separation between two neighboring high- 
lying lines has the asymptotic expression 



Xn-\-l X Tl 



2 

A H n 3 ' 



(146) 



For x > x nov (such that the separation between neigh- 
boring lines becomes small compared to unity), we can 
approximate the sum of optical depths due to Lyman 
transitions as a Riemann integral: 

^2 T rAn{x) ~T C ^2 - x n) -^—e~ {x ~ Xn)2 



- (x ~ u)2 du 



-^erfc(x) =t c - r c (j) c (x). 



Thus, we obtain, using the definition ( |114[ ) 



(147) 
(148) 



Equation ( 121 ) therefore has an analytic solution: 



Vj {x) = ^e-^-^ 2 erfc (x - Xj + ~) ei^+^T 

1 p-(x-Xj) 2 

(149) 

where in the second line we used r c ^> 1 (r c > 10 2 at all 
times). 

We therefore obtain: 



P(j -> i) 



r ± JL e -s(*«-**) 

2 1 



Aid 3 



P(J -»c)»l- 



1 

2J^ 



e 2 



+oo 



2 (Xi Xj ) 



(150) 



erfc(x)e ^ Xi ^dx 



(151) 



where the integral can be evaluated after differentiating 
with respect to xi. The overlap-induced transition rates 
R^Xjp an d photoionization rates can then be ob- 



tained from Eqs. ( p7| and (p9j). 

Similarly, the function v c (x) J (f> c (y)G(x,y)dy satis- 
fies the following differential equation: 



dv c 
dx 



(152) 



which has the solution, valid to lowest order in l/r c : 



v c (x) 



1 

Tc. 



1 — -erfc(x) 



(153) 



20 



The overlap-induced recombination coefficient can then 
be written as: 



(ov) 

OL- = 
W r 2 



8ttz/ c 3 Ah ( 2nh 2 ^ 3/2 



-cr 



TYLp kT\-\~ 



where we used Eqs. (129), (135) and (1381 



(154) 



After some manipulations, we can show that the 
overlap-induced recombination coefficients can be simply 
expressed in terms of the overlap-induced photoioniza- 
tion rates: 



o- (ov) sa 3 
lp V m.fcT, 



2-Kh 2 ^ 3/2 



(155) 



;-3 



From the asymptotic expression of Aj Pt i s oc j u , we can 
also show that 



(156) 



Equations (155) and ([156]) are simply the usual detailed 
balance relations, in the limit vi Vj ss v c . 

The expressions provided in this section are valid for 
n > n ov , when Lyman lines are within less than a Doppler 
width of each other. We therefore expect the expressions 
for p^p* and tti° p v ' to be accurate in the regime where 



they are significant, for n 



> 



^> n ov . On the other 



hand, the interline transition probabilities P(i — > j) 
should smoothly transition from the asymptotic expres- 
sion (150) for i.j > n ov to the Sobolev values for nearly 

i,j < n ov , that is 



instantaneous feedback for 



P(i -> j) = T-^. 



i-r; 



(157) 



We checked that this is indeed the case by integrat- 
ing numerically Eqs. (121) and (125). We therefore set 
P(i j) — f° r min(z, j) < n ov , and use Eq. (150) 



otherwise. Since, as we shall see below, line overlap ap- 
pears to lead to negligible changes to the recombination 
history, the exact value of the interline transition proba- 
bilities near n ov is not critical. 



Results 



As we showed in Section |VD 5| the net rate of change 
of the ground state population through the high-lying 
Lyman transitions vanishes. Since any np — > Is transi- 
tion is systematically followed by the absorption of the 
emitted photon, the high-lying np states are virtually 
radiatively connected to one another (and to the con- 
tinuum), rather than being radiatively connected to the 
ground state. In the language of Ref. [41] , the high- lying 
p states are interior states. The EMLA formalism devel- 
oped in Ref. [UJ can then easily be extended to include 
the overlap-induced transitions, which only depend on 
the matter temperature. 

We added the overlap-induced H(np) O H(n'p) and 
H(np) O e~ + p transitions to our effective rates code. 



We computed the change in the effective recombination 
coefficients A2 S , Ai v and effective 2p —} 2s transition 
rate 7?. 2 p-j-2s when the states 2s, 2p are considered as the 
only interface states (i.e. cutting off the Ly/3 transition 
and above for simplicity; the effect of overlap is indepen- 
dent of that simplification). We find that the relative 
changes in each of the effective coefficients are at most 
a few times 10 -5 . As a comparison, the change in total 
effective recombination coefficient between n max = 100 
and n max = 200 is of order 0.3 to 2% over the tempera- 
ture range considered; from n max = 200 to n max = 400, 
this change is of order 0.08 to 0.6%. Therefore, the effect 
of overlap is a few orders of magnitude smaller than the 
mere error due to the necessary truncation of the high 
energy shells when computing the effective rates. Since 
previous work [521 HO] have shown that MLA computa- 
tions with n max ~ 100 already reach the desired level 
of accuracy, we conclude that line overlap can be safely 
ignored. 



VI. POPULATION INVERSION 

Another important effect to consider is that of popula- 
tion inversion. Sufficiently bottlenecked multi-level sys- 
tems may develop population inversion between radia- 
tively coupled states. If this effect is dramatic, if the ve- 
locity field is sufficiently coherent, and the effective path 
length long enough, stimulated emission may lead to in- 
tense, narrow, coherent maser radiation [68, 69 . Many 
astrophysical masers are now known, and play an impor- 
tant role in firming up the extra-galactic distance scale 
[fjM74] , Most known astrophysical masers are molecular, 
but the emission-line star MWC349 is one example of a 
hydrogen recombination line maser [731 175] . 

It has been suggested both that the recombining pri- 
mordial hydrogen plasma exhibits sufficient population 
inversion for an all-sky natural maser at high n, and 
that rare, extremely overdense regions mase during the 
epoch of recombination [7S1 [77] • The possibility of a 
cosmic recombination maser is particularly enticing, be- 
cause it could conceivably amplify weak low-frequency 
v = 100 Mhz — 1 Ghz CMB spectral distortions from 
recombination to a detectable level [37] . 

Output from RecSparse computations shows that 
some radiatively connected a-transitions (An = ±1) be- 
tween states with I = O (1) and n ~ 50 do show rather 
dramatic population inversion from z < 800 onward. We 
have verified that the width of lines showing population 
inversion is dominated by Doppler and not natural broad- 
ening. If we neglect natural broadening, all pairs of ra- 
diatively connected Al — ±1 pairs for some fixed values 
of n, n' (where n and n' denote the hydrogen shells con- 
nected by the line) have the same line profile, 4> n . n > (x). 
In this case, the steady-state radiative transfer equation 
may be solved to obtain the photon occupation number 
f~ n i on the red side of the line in terms of the occupation 
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number , on the blue side of the liner 



f~ , = f + , + ( f cq , - f + ,) (l - e~ T ».»0 

J n,n •> n,n Wn,n •> n.n J \ } 



T n ,n' — 2_^, T n,n'i 



. l.U ^~n,n f J n,n' 

Jn,n' ^— \ l.l' 



V T 
t-^il.V n.n' 



(158) 
(159) 

(160) 



where T^ l n , is the optical depth in the transition whose 
initial and final states have quantum numbers (n, I), and 
(n' : I'). We use the convention n' > n. The sum is over all 
allowed sets of quantum numbers obeying dipole selection 
rules. The occupation number in equilibrium with this 

transition is /° q „'/ = x n >. v gi/ {x n< ig v ), and gi = 2(21 + 
1) is the statistical degeneracy of an atomic state with 
angular momentum I. 



It is easily seen from Eq. (1581 that significant ampli- 
fication due to stimulated emission requires T n<n i *C — 1. 
We inspect the output of a RecSparse run with 
«max = 180 and find that when r < 0, |r| < 10~ 8 in the 
most amplified case. In the purely radiative case, then, 
there is no high-n cosmological hydrogen recombination 
maser. 

The populations of the 2p and 2s states are strongly 
out of equilibrium towards the end of the epoch of pri- 
mordial recombination, with the ratio X2 P /3x2 S reaching 
values as large as ~ 3 for z ~ 500 — 600. Although 
we do not track the populations of the 2p 3 / 2 and 2p 1 / 2 
states separately, the overlap of the lsi/2 — 2pi/ 2 and 
lsi/2 ^ 2p 3 / 2 Lya doublet ensures that they are in sta- 
tistical equilibrium, so X2 P3/2 /2x2 Sl/2 = x 2p/3x2 S - It is 
therefore a priori possible for the 11GHz 2p 3 / 2 — > 2s 1 / 2 
transition to be amplified. We checked however that 
the largest negative optical depth in this transition is 
t ~ — 10~ 7 for z ~ 1100, and therefore the population 
inversion does not lead to any significant amplification. 



VII. CONCLUSIONS 

In this work, we have evaluated the impact of several 
previously neglected radiative transfer effects on cosmo- 
logical hydrogen recombination: 

• Thomson scattering in the Ly-a line was shown to be 
negligible, with corrections to the recombination history 
\Ax e /x e \ < 3 x 10 -5 . We showed that at early times, 
z > 1300, the dominant effect was a delay of recombi- 
nation, due to the reinjection of photons from the red 



side of the line back into the blue side of the line dur- 
ing large angle scattering events. That effect can only be 
properly accounted for with a full kernel approach, since 
scattered photons are redistributed on frequency scales 
larger than the characteristic scale over which the radi- 
ation field changes. At lower redshifts, recoil becomes 
dominant and Thomson scattering accelerates recombi- 
nation by helping photons escape from the Lya line. 

• Distortions from the deuterium Ly-a line were shown 
to be negligible. Indeed, the very fast D(2p) — > D(ls) 
transition rate, due to the relatively small optical depth 
in the deuterium line, brings the deuterium 2p to Is ra- 
tio close to equilibrium with the incoming radiation field 
on the D Lya line. Moreover, any distortions are further 
washed out due to the very large frequency diffusion rate 
caused by resonant scattering by neutral hydrogen. Ac- 
counting for deuterium therefore leads to changes of at 
most 0(1O~ 5 ) to the recombination history due to the 
small change in the expansion history and the ambiguity 
in defining x e . 

• The high-lying, non-overlapping Lyman transitions 
above Ly7 (strictly) can be artificially cut off without loss 
of accuracy. Only the 2s, 2p, 3p and Ap states therefore 
need to be considered as "interface" states in an EMLA 
computation |41) . Diffusion in Ly/3 and higher lines can 
be neglected. 

• Overlap of the high-lying Lyman lines, as well as over- 
lap of the extremely high-lying lines with the continuum, 
was shown to lead to 0(1O~ 5 ) changes to the effective 
transition rates. This change is a few orders of magni- 
tude smaller than the mere truncation error in the effec- 
tive transition rates computed with n max ~ 100 energy 
shells. The effect of line overlap is therefore negligible. 

• Cosmological hydrogen masers are shown not to arise 
in this purely radiative treatment. 

The goal of the ongoing work in the field is to develop a 
complete theory for hydrogen recombination, with a well 
understood error budget. In this paper, we have eval- 
uated the impact of some radiative transfer effects that 
had not been previously addressed. While it is possible 
that some effects have not been considered yet, we be- 
lieve that most of the radiative transfer effects relevant in 
primordial hydrogen recombination are now well under- 
stood. The picture is less clear for the effect of collisional 
processes [40], for which the rates are relatively poorly 
known. Ultimately we need a recombination code that 
is not only accurate, but also fast, in order to be included 
in Markov chains for cosmological parameter estimation; 
this will be the subject of future work. 
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